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We calculate the unpolarized transverse momentum dependent fragmentation function (TMDFF) 
at next-to-next-to-leading order (NNLO), evaluating separately TMD soft factor and TMD collinear 
correlator. For the first time the cancellation of spurious rapidity divergences in a properly defined in¬ 
dividual TMD beyond the first non-trivial order is shown. This represents a strong check of the given 
TMD definition. We extract the matching coefficient necessary to perform the transverse momen¬ 
tum resummation at next-to-next-to-next-to-leading-logarithmic accuracy. The universal character 
of the soft function, which enters the definition of all (un)polarized TMD distribution/fragmentation 
functions, facilitates the future calculation of all the other TMDs and their coefficients at NNLO, 
pushing forward the accuracy of theoretical predictions for the current and next generation of high 
energy colliders. 


Introduction. Multi-differential cross sections play 
a central role in our understanding of QCD dynamics. 
In this context the definition of transverse momentum 
dependent functions (TMDs) has been recently revis¬ 
ited, updating the pioneering work of Collins and Soper 
[1, 2], in order to solve the subtle issue of the cancella¬ 
tion of spurious rapidity divergences inside an individual 
TMD. As a result, one has achieved the factorization 
theorems for Drell-Yan, Vector Boson/Higgs Produc¬ 
tion, Semi-Inclusive Deep Inelastic Scattering (SIDIS) 
and e+e“ -T 2 hadrons processes in terms of individ¬ 
ually well-defined TMDs [3-5]. All these processes are 
fundamental for current high energy colliders, like the 
LHC, KEK, SLAC, JLab or RHIC, and future planned 
facilities, like the EIC, AFTER, the LHeC or the ILC. 

While the formulation of the factorization theorems is 
solid, a direct evaluation of an individual TMD at two 
loops is still lacking. Such a calculation provides a funda¬ 
mental check of the factorization theorem and important 
information for data analysis. The one-loop TMDs with 
various quantum numbers have been computed by sev¬ 
eral groups [5-11]. At two loops some properties of TMDs 
have been deduced from cross section calculations carried 
out in QCD (see e.g. [12-16]). In this work we present 
the result of the next-to-next-to-leading-order (NNLO) 
evaluation of the unpolarized transverse momentum de¬ 
pendent fragmentation function (TMDFF), which is an 
essential part of TMD factorization theorems, both in the 
non-singlet and singlet channels. 

The evaluation of individual TMDs at higher orders 
is an utterly non-trivial check for their definition, since 


starting from the two-loop order, the singularities of var¬ 
ious types mix up. The two-loop calculation (for the first 
time presented in this letter) shows that the combination 
of factors disentangles and cancels the spurious rapidity 
divergences within the proper definition of one TMD. The 
intermediate pieces of the calculation are also relevant 
per se. In fact we define a regulator for rapidity diver¬ 
gences which can be used in combination with standard 
dimensional regularization for evaluation of any TMDs. 
The Soft Function (which is essential for the TMD def¬ 
inition and whose two-loop result will be presented in 
a forthcoming publication [17]) is a key element for the 
NNLO calculation of all (polarized) TMDs. 

In this letter we present the matching coefficient of the 
unpolarized quark TMDFF onto the integrated fragmen¬ 
tation function (FF) at NNLO (both non-singlet and sin¬ 
glet channels), using explicitly the formalism of Ref. [5]. 
This result can be immediately used in forthcoming phe¬ 
nomenological applications, see e.g. Ref. [18, 19]. Our 
consideration fills the gap in TMD phenomenology, be¬ 
cause NNLO coefficients for Transverse Momentum De¬ 
pendent Barton Distribution Function (TMDPDF), an¬ 
other important ingredient of TMD factorization theo¬ 
rems, can be extracted from the NNLO calculations made 
in the related formalism [12-16]. 

Definitions. Individually TMDs are defined as a 
product of two separate matrix elements, the (square 
root of the) Soft Function (SF) and the corresponding 
collinear matrix element. The SF is a spin and process 
independent vacuum expectation value of Wilson lines 
[3, 4], and it is defined as 
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where S'„ and Sn stand for soft Wilson lines along light-cone directions^, and n and n are light-cone vectors {rP = 
rl^ = 0, n ■ n = 2). The unpolarized collinear matrix element is defined as 
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where index i refers to the parton flavor, Ph is the hadron 
momentum, z is the Bjorken variable, = P^ jz and 
fc„_L = —Ph^jz. The subscript “zb” stands for zero-bin 
subtracted^. See Ref. [5] for more details regarding the 
particular definition of Wilson lines. 

Individually both matrix elements have rapidity di¬ 
vergences at every order in the perturbative expansion. 
These divergences are neither ultraviolet (UV) nor long¬ 
distance ones and, in principle, are not sensitive to con¬ 
fining dynamics [3, 4, 7, 27]. As argued in Refs. [3, 4, 7], 
such divergences can be removed in the correct combina¬ 
tion of soft and collinear matrix element. 

The essential property of the SF which allows one to 
remove the rapidity divergences, is that the logarithm of 
the SF is maximally linear in the logarithmical rapidity 
divergences. Therefore, it can be split into two pieces [4], 

5(A, (A, L5+/,) (A, L,s -), (3) 

where v is an arbitrary, dimensionless and positive real 
number that transforms as p'^ under boosts and ^ we 
introduce the convenient notation 

hx = ln(AA2e2T'®/4). 

Variables are rapidity regulators that one uses in the 
n- and n-collinear sectors (our implementation of it is 


^ The superscript T on Wilson lines in Eq. (1) implies subsidiary 
transverse links from the light-cone infinities to transverse in¬ 
finity, see details in Refs. [20-22]. These links guaranty gauge 
invariance and are necessary for calculations in singular gauges. 
The presented calculation has been performed in Feynman gauge, 
where the contribution of transverse links vanishes. 

^ The zero-bin subtraction is the term used in the Soft Collinear 
Effective Theory (SCET) literature to account for the double 
counting with the soft sector. These definitions are equivalently 
stated in QCD and SCET, see e.g. Refs. [4, 23—26]. Here, we ap¬ 
ply this term since the definition that we follow [5] was originated 
within SCET. However, the present calculation is performed in 
standard QCD. For the used regularization, the application of 
zero-bin subtraction is equivalent to calculate the colliner ma¬ 
trix element naively (Eq. (2)) and then subtract the soft function 
matrix element in Eq. (1), thereby obtaining the so-called “pure 
collinear” matrix element: Apure ~ Anaive5'“^. The precise de¬ 
tails on the definition will be presented in Ref. [36] . 

^ We denote by and p~ the large components of the incom¬ 
ing and outgoing parton momenta, respectively, in a SIDIS hard 
process. 


specified later in Eq. (6)). Tildes mark quantities cal¬ 
culated in the coordinate space. In our calculation the 
relation in Eq. (3) has been checked explicitly at NNLO. 

The result of the combination of one piece of the SE 
and the collinear correlator (A) is free from rapidity di¬ 
vergences and hence can be considered as a valid hadronic 
quantity. Eor the unpolarized TMDEE in coordinate 
space we have 

Di-^h{z, J-ifi, Ii^q) = Aj^^(z, L^, ) 55 (L^, L,5+/j^) 

= ^i^h{z,'Lfj^,Xg-) S (]4) 

where we have introduced the shorthand notation 

lx = ln(A/V), \s- = ln{6~/p~). 

In this equation Ai^h represents the naively calculated 
collinear matrix element, with no subtraction of the over¬ 
lapping with the soft region. If the hard scale in the pro¬ 
cess, say the mass of the virtual photon in e+(p)e“(p) —>■ 
2 hadrons, is given by Q^, then A Cd are fractions 
of Q^, satisfying (fCd = where A = and 

A = (ia the following we omit the subscripts 

F, D where unnecessary). At small values of the impact 
parameter b the renormalized TMDEE can be factorized 
again in 

~ dr ! z \ 

L/ij Ic) “ y j-3—2e j 

xdj^h{T,p), (5) 

where di^hi^-, p) is the renormalized integrated FF. In 
Eq. (5) and in the rest of this letter, the repeating flavor 
index implies summation. The outcome of this work is 
the calculation at NNLO of the non-singlet and singlet 
part of quark to quark and quark to anti-quark coeffi¬ 
cients, respectively Cq^q{z, L^, 1^), Cq^^q^ {z, L^, 1^) and 

Cq—iq{z, L^, 1(^). 

Regularization. The choice of the infrared (IR) and 
rapidity regularization scheme is one of the central points 
for the evaluation of TMDs. The regularization should 
satisfy several important demands, such as: it should re¬ 
spect the exponentiation property of Wilson lines (which 
is necessary for the existence of the relation in Eq. (3)); 
it should match the singularities of the naively calcu¬ 
lated collinear matrix element in the soft limit with the 
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ones of the SF (which is necessary for a straightforward 
treatment of the zero-bin subtraction, see footnote 2). 
Additionally, the chosen regularization scheme should be 
convenient for multi-loop integral computations. 

One of the popular choices of regularization is to use 
tilted Wilson lines, see e.g. Refs [6, 9]. However, with 
this regularization the number of loop-integrals and their 
difficulty is significantly higher than with others. The 
analytical regulator, that was used in NNLO calculation 
in Refs. [15, 16], is highly efficient for computation and 
satisfies the necessary requirements, but it is not capable 
of regularizing rapidity divergences of the SF, which is 
crucial in the proper definition of an individual TMD. 

Here we regularize the rapidity divergences with the S- 
regularization, that has been used for the same purpose 
by many authors, see e.g. Refs. [5, 7, 8, 28]. To regularize 
the rest of UV and IR divergences we use standard di¬ 
mensional regularization (DR) with D = 4 —2e, while the 
incoming/outgoing partons are on-shell and massless 
To match the required demands at multi-loop level the 
(5-regularization is here modified. First, in order to sup¬ 
ply the non-abelian exponentiation property [29, 30], and 
hence the relation in Eq. (3), the (5-regulator should be 
implemented at the operator level, see e.g. the discussion 
in Ref. [31]. We thus modify the definition of Wilson lines 
as ^ 

Wfi{0) = Pexp 
—>■ P exp 
Sn{0) = Pexp 
—Pexp 
S'„(0) = Pexp 
—>■ P exp 

where (5^ —?► 0+. Second, in order to match the IR soft 
singularities of the naively calculated collinear matrix el- 
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^ For renormalization we use MS-scheme with the rescaling factor 

(47re^E)^ 

^ We should mention that the presented regulator has some in¬ 
conveniences typical of such regularizations. One of them is the 
potential violation of gauge invariance. However, artificial gauge 
violating terms can be easily traced and discarded. Another in¬ 
convenience is that 5 parameter regularizes not only rapidity di¬ 
vergences but also some other soft divergences. In general, this is 
not a problem, since all soft divergences cancel in the final result. 
Nonetheless, we made a complete analytical calculation, where 
different soft divergences have different signature, and checked 
the cancellation individually for every sector. The details will be 
presented in Ref. [36]. 


ement and the SF, the S in Ai^h should be rescaled by 
2 , i.e. (5 —>■ 6/z. 

Such modified regularization is appropriate for being 
used in multi-loop calculations and for the evaluation of 
the relevant matrix elements separately. 

Extraction of the matching coefficient. In order 
to extract the matching coefficient at NNLO one needs 
to evaluate the SF and the collinear matrix element in 
Eq. (2) partonically at NLO and NNLO. The obtained 
functions, together with the renormalization multipliers 
should be combined into the partonic TMDFF, Eq. (4). 
The partonic TMDFF then is matched onto the inte¬ 
grated FF in the operator product expansion sense. 

At one-loop order the procedure is presented e.g. in 
Ref. [5] . The complete expression for TMDFF reads 

Dfi, = A^X - ^ - 4 " + 4 '> ( 7 ) 

where Z 2 is the quark wave-function renormalization 
constant, and Zd is the TMDFF operator renormaliza¬ 
tion constant. Note that in this expression, as well as 
in Eq. (10), zero-bin subtractions are explicitely taken 
into account (the SF is subtracted instead of added). 
Throughout the paper we use superscripts in square 
brackets to denote the order in the perturbative expan¬ 
sion, e.g. S = where Os = j£pj^ and also the 

shorthand z = 1 — z. 

The rapidity divergences appear in both and 
but cancel in Eq. (7). The ultraviolet divergences are 
renormalized by Z 2 and the suitably chosen Zjj. There¬ 
fore D^q^q is a function of z, L^, f^ and e, which regu¬ 
larizes the IR collinear divergences. The collinear diver¬ 
gences are part of the integrated FF, while the matching 
coefficient C is given by 


rW = nW - 
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At one-loop order we obtain the well known result [3, 5] 
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and the trivial result = 0. Here, Vq^q{z) = 

((1 -I- z^)/z)+ is the quark splitting function. The plus- 
distribution is defined as {f{z))^ = f{z) — 5{z) dyf{y). 
At two-loop level the TMDFF is 
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The two-loop rapidity divergences appear only in the first 
line of Eq. (10). Notice that in contrast to NLO, where all 
rapidity divergences arise with S{z) prefactor and cancel 
trivially between SF and A, at NNLO the rapidity di¬ 
vergences arise with an involved 2 -dependent structure. 
At two-loop level the rapidity divergences of and 
^[21 

mix up with UV divergences, and the mixture can¬ 
cels in the combination in Eq. (10). In general, Eq. (10) 
possesses a complex system of cancellations of various 
divergences [36]. The realization of all this cancellation 
represents an important check of our calculation. 

The matching coefficient at two-loop level is given by 
the combination 


Z—>■/ l—¥f i—¥k 
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where the symbol ® denotes Mellin convolution in the 
Bjorken variable z, while A: is a flavor index. Clearly, 
each addend of this sum is free of rapidity divergences. 

Renormalization group and matching. The 
renormalization group equations of the TMDFF and the 
integrated FF provide also important checks for our cal¬ 
culation. We have that 

ID = K^spk - ik , (12) 


where Tcusp is the cusp anomalous dimension, e.g. 
Ref. [4]. The result for jy is extracted from the calcula¬ 
tion of the non-singlet part of the quark form factor [32] . 
Then we have 

= (13) 

which allows the resummation of the rapidity logarithms. 
Putting together Eqs. (5), (12) and (13), one finds 
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where the convolution is understood in the Bjorken vari¬ 
able 2 and 

ICU,{z) = - ikm) - (15) 

The function 'Pi^j{z) is the DGLAP kernel for the inte¬ 
grated FF at NNLO (see e.g. Refs. [33-35]). 

The evolution equations allow one to write the Wilson 
coefficient in a more compact form; 

= exp [ - P*L^] . (16) 

The most general structure of Cy of the n-th perturbative 
order is 


2n 




(n;fe) 

ij 


T ^ 

■ 


k=0 


(17) 


The coefficients are related by the recursive rela¬ 

tion 
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Thus, given the expressions for the anomalous dimensions 
one needs only the boundary coefficients in order 

to reproduce the complete expression for the matching 
coefficient. In our calculation we evaluate the complete 
logarithmical structure of the TMDFF and explicitly con¬ 
firm the relations in Eqs. (16) and (18), thus providing a 
strong check for the whole calculation. 

Results. For completeness we present LO and NLO 
expressions for boundary conditions. They are 


A(0;0) 

rim 
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^ (4p(2)lnz -b 22 )+ -b (5(2) ^6 - , 


where p{z) = The corresponding quark-antiquark 

coefficients are zero. 

The NNLO coefficient can be decomposed as 
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Then the functions Qi are 
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For the 2-loop singlet part we obtain 
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These expressions represent the main result of this letter. 

Conclusions. TMDs are defined as the product in 
coordinate space of the collinear matrix element and the 
square root of the soft function. In this paper we pro¬ 
vide the explicit check of this statement for the first time 
at NNLO, for the quark TMDFF. The calculation (per¬ 
formed within standard QCD and in Feynman gauge) in¬ 
cludes the independent computation of the soft function 
and the collinear matrix element, and their subsequent 
recombination into a well-defined TMD. We have refor¬ 
mulated the IR and rapidity regularization of Ref. [5] in 
order to extend the definition of an individual TMD to 
multi-loop level. We obtain the complete analytical ex¬ 
pression for the TMDFF, and comprehensively investi¬ 


gate the structure of soft/rapidity singularities and their 
cancellation. The cancellation of singularities provides 
a strong check of the final result. As a further check 
we find a complete agreement between the logarithmi- 
cal part of the final result and the known predictions 
of renormalization group. The soft factor that has been 
evaluated in this work, is universal and spin-independent, 
and thus can be used for the calculation of all TMDs at 
NNLO. Finally, the calculation of the TMDFF performed 
in this work allows us to extract the relevant perturba¬ 
tive matching coefficient at NNLO, necessary to perform 
the resummation of large logarithms at NNNLL, pushing 
the phenomenology a step forward. The applied method 
can be readily used to obtain other relevant perturbative 
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ingredients. The detailed report, including the other fla¬ 
vor parton contributions, the explicit expressions, as well 
as the description of the calculation, will be given in a 
separate publication [36]. 
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Changes in arXiv version 3 compared with version 2 

Equation (21) has been corrected, and it is consistent with Ref. [36]. 



